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Electromagnetic bound states in the radiation continuum are studied for periodic double arrays 
of subwavelength dielectric cylinders in TM polarization. They are similar to localized waveguide 
mode solutions of Maxwell's equations for metal cavities or defects of photonic crystals, but, in 
contrast to the latter, their spectrum lies in the radiation continuum. The phenomenon is identical 
to the existence of bound sates in the radiation continuum in quantum mechanics, discovered by von 
Neumann and Wigner. In the formal scattering theory, these states appear as resonances with the 
vanishing width. For the system studied, the bound states are shown to exist at specific distances 
between the arrays in the spectral region where one or two diffraction channels are open. Analytic 
solutions are obtained for all bound states (below the radiation continuum and in it) in the limit of 
thin cylinders (the cylinder radius is much smaller than the wavelength). The existence of bound 
states is also established in the spectral region where three and more diffraction channels are open, 
provided the dielectric constant and radius of the cylinders are fine-tuned. The near field and 
scattering resonances of the structure are investigated when the distance between the arrays varies 
in a neighborhood of its critical values at which the bound states are formed. In particular, it is 
*^ ■ shown that the near field in the scattering process becomes significantly amplified in specific regions 

of the array as the distance approaches its critical values. The effect may be used to control optical 
non-linear effects by varying the distance between the arrays near its critical values. 
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I. INTRODUCTION 

Scattering problems in quantum mechanics and Maxwell's theory exhibit a great deal of similarities. The dielectric 
constant of an electromagnetic scattering structure is analogous to a scattering potential in the Schrodinger equation. 
If the potential is real and bounded from above, then solutions of the wave equation (Schrodinger's or Maxwell's) 
appear in two kinds [2|. There are radiation modes that carry the energy flux to the spatial infinity and modes 
trapped in some (bounded) spatial region. The former are known as scattering states. They do not have a finite L2 
norm, and their spectrum is continuous (the radiation continuum). In contrast, solutions of the second kind have a 
finite L2 norm and a discrete spectrum which typically lies below the radiation continuum. In quantum mechanics, 
these are bound states. In Maxwell's theory, they are electromagnetic waves localized, e.g., in metal cavities or in 
defects of photonic crystals. It was first proved by von Neumann and Wigner [l| in 1929 that, in quantum mechanics, 
there may exist bound states in the radiation continuum (i.e., there are localized solutions of the Schrodinger equation 
whose energies lie in the continuum part of the system spectrum). An apparent physical peculiarity of such states 
is that, despite the boundedness of the potential barrier, they do not decay into the scattering modes through the 
tunneling, which is rather counterintuitive and unusual. An example constructed by von Neumann and Wigner was 
somewhat artificial from the physical point of view, and no significant applications of this phenomenon immediately 
followed their discovery. However, much later unusually long-lived resonances were observed in atomic physics [3| and 
the question about a physical mechanism of their stability arose. 

From the mathematical point of view, scattering resonances can be associated with solutions of the wave equation 
that satisfy the outgoing wave boundary condition at the spatial infinity, i.e., in the asymptotic region, the solu- 
tions look like waves outgoing from the scattering structure 0j. Clearly, such solutions cannot preserve energy (or 
probability) flux and, hence, do not exist for real frequencies (or energies). So they have complex eigen- frequencies 
and decay exponentially with time. Due to linearity of the wave equation, a solution of a scattering problem is a 
linear superposition of the incident wave and scattered waves that satisfy the outgoing wave boundary condition. The 
resonant solutions may therefore be used to construct the scattered wave, and the flux conservation is achieved by 
balancing the incoming fiux of the incident wave and the outgoing flux carried by scattered waves. It was proved 
then that the resonant scattering amplitude may be approximated by the celebrated Breit- Wigner formula (2i 14] 
A(uj) ~ iT/[cu — ujr + *r) where uj is the incident wave frequency, and w,- — iV is the eigen- frequency of the resonant 
mode. If the imaginary part F is small {cUr ^ F), then the scattering cross-section a{uj) ^ \A{uj)\'^ exhibits a sharp 
maximum (of width F) at frequencies near w^. Thus, resonances with small F can be associated with quasi-bound 
states, excitations with a finite lifetime ~ 1/F, that decay by emitting a nearly monochromatic radiation of frequency 
LOr- In this picture, the bound states in the radiation continuum may be understood as scattering resonances with the 
vanishing width. If a system allows for two close resonances whose positions and widths depend continuously on a 
physical parameter of the system that can be varied, then there might exist critical values of this parameter at which 
one of the resonances turns into a bound states p-01 ■ This is the subject of the theory of coupled resonances developed 
in quantum scattering theory [^ . The conditions under which the width of one of the coupled resonances vanishes 
were established. The theory of coupled resonances provided an explanation of the existence of the aforementioned 
unusually long-lived resonances observed in atomic systems. 

Resonant scattering properties of gratings and periodically structured thin films were observed long time ago in 
optics [8l-[l0l|. Typically, a grating would exhibit sharp resonances near its diffraction thresholds [IJ. However an 
explanation based on quasi-bound states came much later. Naturally, the Breit- Wigner theory can be extended to 
the electromagnetic scattering, and coupled resonances might be expected to exist in photonic structures jl2,,, 17]. A 
simple example is provided by two parallel identical gratings. Each grating supports a quasi-stationary electromag- 
netic mode localized in its vicinity (a trapped mode) , which can be excited by an incident radiation. In the zero order 
diffraction, this mode has only one decay channel. Therefore one might anticipate a formation of a wave guiding 
mode confined between the structures at a certain distance between the gratings. If the trapped mode associated 
with each grating decays by emitting a monochromatic radiation outgoing both ways along the normal to the struc- 
ture, then there might exist a distance at which the waves outgoing from each structure interfere destructively in the 
asymptotic region, thus producing a localized stable solution (the energy flux outgoing from the structure vanishes). 
This argument can be made rigorous by considering the Fabri-Perot interferometer whose interfaces have reflection 
and transmission coefficients described by the Breit- Wigner theory [Ij, [l3]. When the distance gets smaller, the 
interference effects of the near field of the trapped modes becomes significant, and the Fabri-Perot argument becomes 
inapplicable as it applies only to the modes carrying the energy fiux, i.e., radiation modes. Maxwell's equations must 
then be studied. Bound states in the radiation continuum were indeed found by numerical studies of the system of 
two gratings in the spectral region where only one diffraction channel is open. They are formed at a discrete set 
of distances that follows the pattern predicted by the aforementioned Fabri-Perot argument for large distances and 
substantially deviates from it when the distance is of order the grating period [12j . Also, localized waveguide solutions 
were argued to exist in the same spectral region by studying Maxwell's equations for layers whose dielectric function 



is periodic in one direction along the layer and translational invariant along the other one |18| . 

Here bound states in the radiation continuum are studied analytically for a system of two arrays of parallel dielectric 
cylinders. The approach is based on the resonant scattering theory [2, [l3] where the bound states are identified as 
resonances with the vanishing width (the distance between the arrays is a physical parameter which regulates the 
coupling of the resonances). The study is carried out for the whole radiation spectrum range. In addition to the 
bound states in the zero-order diffraction (which were known to exist from the early numerical studies) , bounds states 
are found in the spectral range where two diffraction channels are open. Analytic solutions of Maxwell's equations for 
all the bound states are given in the limit when the cylinder radius is much smaller than the period of the structure. 
The system is shown to have bound states below the radiation continuum whose explicit form is also found. 

The interest to subwavelength periodic structures has been revitalized due to a recent technological progress in man- 
ufacturing such structures for the optical frequency range [ij]. An ultimate goal is to develop a technology for the 
so-called all optical data processing [i5l[i6|. One way or the other, this task amounts to finding physical mechanisms 
to control optical non-linear effects in photonic structures as, e.g., a typical logical element in data processing requires 
a device similar to a transistor. As a potential application, here the near field of resonant excitations is studied as a 
function of the distance between the arrays that varies in the vicinity of its critical values at which the bound states 
are formed (the width of one of the resonances vanishes) . It is shown that the near field can be significantly amplified 
(as compared with the incident radiation field) in some particular regions of the structure as the distance approaches 
its critical values. The effect can therefore be used to enhance optical non-linear effects in a controllable way. 

The paper is organized as follows: a classification of solutions of Maxwell's equations that is convenient for identify- 
ing the bound states in the radiation continuum in the system studied is described in Section |TT1 Explicit solutions for 
the bound states below the radiation continuum as well as in the radiation continuum when one diffraction channel 
is open are given in Section IIIII The effect of the near field amplification for a resonant scattering when the distance 
between the arrays is close to one of its critical values at which the bound states are formed is investigated in Sec- 
tion nil CI Lastly, the existence of bound states when more than one diffraction channel are open is established in 
Section |IVl Whenever possible, the technical details are moved to the appendix. 



II. SCATTERING THEORY AND CLASSIFICATION OF THE FIELDS 

The system considered is sketched in Fig. [TJa). It consists of an infinite double array of parallel, periodically 
positioned, dielectric cylinders suspended in the vacuum [1^ . The cylinders are assumed to be non-dispersive with a 
dielectric constant Sc > 1. The coordinate system is set so that the cylinders are parallel to the y— axis, the structure 
is periodic along the x— axis, and the z— axis is normal to the structure. The unit of length is taken to be the array 
period and the arrays have a relative mismatch a G [0, i] along the x-axis. 
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FIG. 1. Panel (a); Double array of dielectric cylinders. The unit of length is the array period. The axis of each cylinder is 

parallel to the y-axis and is at a distance h from the x-axis. 

Panel (b): The energy spectrum for a Schrodinger equation with radially symmetric potential consists of a discrete spectrum 

of bound states with negative energies and the radiation continuum. The latter may contain additional bound states; these are 

the bound states in the radiation continuum. 

Panels (c) and (d): Harmonically time dependent solutions to Eq.([l| are labeled by points of the spectral cylinder as explained 

in text. The diffraction thresholds £±n are rings that partition the cyhnder into sections corresponding to a fixed number of 

open channels. In particular, £o is the threshold for the radiation continuum below which there can be no scattering states. 



The structure is illuminated by a linearly polarized monochromatic beam with the electric field parallel to the 
cylinders (TM polarization). In the given settings, Maxwell's equations are reduced to the scalar wave equation for a 
single component of the electric field, denoted E, in the x, z— plane: 

^d^E -AE = 0, (1) 

where the dielectric function e is piecewise constant; it is equal to Ec > 1 on the scatterers and to 1 otherwise. For a 
harmonic time dependent electric field E{r,t) = E^{v)e~^^* , the amplitude -Ew(r) satisfies the equation 

HE^ = [-A + k\l-e)]E^ = k'E^, (2) 

where k^ = uP' j (? is the spectral parameter. This equation is similar to the Schrodinger equation with an attractive 
(negative) potential 1 — £. If k is the wave vector of the incident radiation, then the periodicity of the structure 
requires that solutions to Eq.® satisfy Bloch's theorem, 

£;^(r + ei) = e^'^^S^(r), (3) 

where k = fc^-ei + k^e.^^ and {ei, 62, 63} are the unit vectors along the coordinate axes. Consequently, the wave vectors 
of the scattered radiation can only have the following form 

k± = (fc, + 2^m)ei ± Vfc2 - (fcx + 2^TO)2e3 , m = ... - 2, -1, 0, 1, 2... (4) 



with the convention -^/fc^ — (k^ + 27rm)2 = i^ (kx + 27r?7i)2 — fc^ if fc^ < (k^ + 27rm)2. The values of m for which 
^^ > (fca; + 27rTO)2 label a'pen diffraction channels, all other values of m correspond to closed channels. If, for instance, 
the incident radiation propagates in the positive z— direction as shown in Fig. [TJa), then the transmitted field reads 

£;„(r) = e*'^•"■ + ^^,„e*'^-•^ z>h + R, (5a) 

rn 

where T^ are the amplitudes of the transmitted modes and the first term describes the incident radiation whose 
amplitude is set to one. Similarly, in the region to the left of the structure, the electric field is given by 

^^(r)=e*-^ + ^i?,„e"'™% z<-h-R, (5b) 

rn 

where Rm are the amplitudes of the reflected modes. The choice between k+ and k^ is dictated by the outgoing wave 
boundary condition at the spatial infinity |2;| -> 00. In particular, the field corresponding to closed channels in Eqs.® 
decays exponentially in the asymptotic region |z| — > 00. So the closed channels do not contribute to the energy flux 
carried by the scattered wave. The scattered flux is carried only by the radiation in open diffraction channels. The 
condition fc^ = {kx + 27rm)2 defines the threshold for opening the to*'' diffraction channel. The coefficients Rm and 
Tm for open channels are the reflection and transmission coefficients, respectively. 

The transmission and reflection amplitudes may be inferred from the solution _E^ of the Lippmann-Schwinger 
integral equation [2] 

E^r) = e'''"- + ^J (£(ro) - 1) EUro)G{r\ro)dr, , (6) 

in which G(r|ro) = iTTHg{k\r — r(,|) is the 2-dimensional free-space Green's function for the Helmholtz operator A + k'^ 
with outgoing boundary conditions, and Hg is the Hankel function of the first kind of order 0. It should be stressed 
that this equation is to be understood in the distributional sense. This is because the potential under consideration 
is neither compactly supported nor does it vanish at infinity so that the usual methods that establish this equation 
cannot be used ^2||. See Appendix 1X1 for details. 

On a different note, observe that, unless there is no incident radiation (the term e^^'' is omitted), the solutions 
to Eq.(l2]) cannot be square integrable along the z-axis. By the analogy with quantum scattering theory for radially 
symmetric potentials, such solutions are in the radiation continuum of the energy spectrum. In quantum theory, 
solutions to the Schrodinger equation are eigenvectors of the energy operator whose spectrum £ form a subset in 
the real line. The upper part (positive) of the spectrum is continuous and corresponds to scattered states that can 
carry the probability flux to the spatial infinity. Below the continuous part of the spectrum, i.e., £ < 0, the energy 
spectrum is discrete (see Fig.[Hb)). It corresponds to bound states. Bound states have a finite L2 norm (they decay 
fast enough at the spatial infinity). It was first proved by Von Neumann and Wigner in 1929 that under special 



circumstances, there might exist bound states in the radiation continuum. A counterintuitive physical pecuharity of 
such states is that a bound state is a standing wave in a potential well (an attractive potential), while the conventional 
quantum mechanical wisdom would suggest that for a potential bounded above a standing wave with the energy in 
the radiation continuum should tunnel through the potential barrier to the spatial infinity and, hence, cannot be 
stable. Nevertheless, such states do exist and the theory of their formation is now well developed [J,[a, 0|- The goal 
here is to establish a similar picture for electromagnetic excitations in the periodic double array of dielectric cylinders, 
and, specifically, to find the conditions on the physical parameters of this system under which the bound states in 
the radiation continuum exist, their eigen-frequencies, and the analytic form of the corresponding electromagnetic 
fields. The very existence of bound states for this system was first demonstrated by numerical simulations [12]. Here 
a complete analytic study of the system is given. 

Any bound state is a solution of the generalized eigenvalue problem ^ (no incident radiation term) and, hence, is 
fully characterized by the spectral parameter £ — k'^ > and the Bloch phase factor e'^''" because of the boundary 
condition ([3]). The pair (f ,e*'^=") is viewed as a point on the (half) cylinder R+ x S^ where £ € K+ and e^'^" £ S^. 
It will be called a spectral cylinder (or spectral space) of a periodic grating. In contrast to quantum mechanical 
bound states in spherically symmetric systems, the Bloch boundary condition requires a more adequate classification 
of bound states here. In order to identify bound states in the radiation continuum, one has first to determine the 
region of the spectral cylinder occupied by the radiation states. In the asymptotic region \z\ — > cx), a harmonic time 
dependent solution E = Ei^e~^'^* to Eq.(IT]) is characterized by the pair (£, fc^) where £ — k^ — txP'lc?. Given £ 
and k^^ the field outside the scattering region is completely determined by its behavior in the diffraction channels 
as specified in Eqs.®. Consequently, if two radiation modes have the same f , while their parameters kx differ by 
a 27r-multiple, say, by 27rTOo, then they have exactly the same open diffraction channels because in the classification 
introduced in Eq. (j?]) the difference of channels would merely mean the relabeling m ^ m — mo (or the same change 
of the summation index in Eq. (O). Therefore the radiation modes correspond to points (f ,e*'^'=") on the spectral 
cylinder for which one or more diffraction channels are open. 

The spectral cylinder can be partitioned into sections associated with a fixed number of open diffraction channels. 
The diffraction thresholds appear as curves separating these portions of the cylinder. Indeed, let [k^] designate the 
argument of e**^^ in (— tt, tt], i.e., [kx] = kx mod 2tt. Then, up to the aforementioned reordering, the diffraction 
thresholds on the spectral cylinder are exactly 



and they appear in the order. 



£±n{kx) = {2TTn±\[kx]\)\ n = 0,l,2,3..., (7a) 



£o < £-1 <£i< £-2 <£2<£-3<£3--- (7b) 



If kx is identified with the angular variable spanning the compactified direction of the spectral cylinder, then the 
diffraction thresholds are curves in an ever rising order on the spectral cylinder with nodes on the axes [kx] = and 
[kx] = TT. The curve £ = £o{kx) is the threshold below which no radiation modes exist, and therefore, the radiation 
continuum lies immediately above this curve. This continuum is split into distinct regions corresponding to a fixed 
number of open channels by consecutive thresholds as indicated by Eq. (j7bp . These regions will be labeled as radiation 
continuum I, radiation continuum II, radiation continuum III, etc., where the Roman numeral indicates the number 
of open channels in each region. See Fig. [TJb), (c), and (d). 

All the solutions to Eq.Q below the threshold £o, if any, must be bound states, i.e., they decay exponentially in 
all diffraction channels and, hence, have a finite L2 norm in the space 5*^ x R spanned by (x, z) {x is compactified 
into a circle S^ because of the boundary condition ([3])). In contrast, radiation modes behave as harmonic functions 
in the asymptotic region \z\ — >■ cx) and, hence, do not have a finite L2 norm. So the problem is to find, if any, 
square integrable solutions on S*^ x R above the curve £ = £0 which are the sought-for bound states in the radiation 
continuum. 

III. BOUND STATES 

As defined above, bound states are square integrable solutions of the homogeneous Lippmann-Schwinger integral 
equation 

E^{t) ^ — {e{T,) ~ I) E^{v,)G{T\T,)dv, , k^>0, (8) 

which satisfy Bloch's boundary condition ([3]). The square integrability here means a finite L2 norm in 5^ x M (with 
the x-direction compactified into a circle). As also noted above, Eq.([7]) is understood in the distributional sense (see 



Appendix 1X1 for details). This is a generalized eigenvalue problem because the Green's function G'(r|ro) also depends 
on the spectral parameter k^. In this problem, the Bloch parameter k^ may be restricted to the interval [— 7r,7r]. It 
must be stressed that this restriction is only permitted by the absence of an incident wave which otherwise determines 
the phase factor in ([3]). 

The task is to determine the values of a, k^, h and k that allow for the existence of nontrivial solutions E,^ to 
Eq.® for fixed radius R and fixed dielectric constant Sc in the limit of thin cylinders, i.e., kR <C 1. In what follows, 
the existence of a bound state always means the existence of a wavenumber k at which the bound state E^^ occurs. 
The main results established in the present study may be split into the following cases which are also summarized in 
Fig. Ha): 

• Below the radiation continuum: Bound states exist for all k^, a and for all distances between the arrays. 

• Continuum I (one open diffraction channel): Bound states only exist if 

— either kx = and a G [0, i] is arbitrary 

— or a G {0, i} and k^ G (— "■, ■"') is arbitrary 

Under these conditions, for each pair {a,kx) there is a discrete set of distances between the arrays at which 
bound states exist. 

• Continuum II (two open diffraction channels): Bound states only exist for 
a — Q OT a — \ and for a certain dense set of values of kx ■ 

For each allowable pair (a, fc^:), there is exactly one or two distances between the arrays at which the bound 
states exist. 

• Continuum N, N> 3 (three or more open diffraction channels) : Bound states exist only for specific values of the 
radius R and the dielectric constant En- 
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FIG. 2. Panel (a): Values of the parameters a and kx for which bound states exist. The bottom level represents the modes below 
the radiation continuum, the second level represents the modes in the radiation continuum I, and the third level represents the 
modes in the radiation continuum II. 
Panels (b) and (c): The shaded areas represents the region of validity of inequalities (|29|l and (|39|l . respectively. 



In all the above cases, bound states occur in two types, symmetric and skew-symmetric relative to the transformation 
z — ^ —z and x — ?> a+x, when the shift parameter takes the boundary values, a = or a = i. However, this classification 
cannot be established for double arrays with intermediate values of the shift parameters, a G (0, i). 

The above classification of the bound states in the cases a = and a = i is relevant when analyzing the (discrete) 
values of the distance 2h between the two arrays at which the bound states occur. When the double array is symmetric 
(a = 0), then these values always exist in the whole range of /i G (i?,oo). On the contrary, skew-symmetric bound 
states will only occur for those values of h that exceeds a certain minimal threshold, and this threshold increases as 
R^{ec — 1) — > 0. In other words, when the scattering structure becomes more transparent, the two arrays have to 
be taken further apart in order for skew-symmetric bound states to form. A similar phenomenon is observed for the 
skew-symmetric array (a = i). In this case, however, the skew-symmetric modes behave as the symmetric ones in 
the previous case and vice versa, i.e., the minimal threshold distance exists for the symmetric modes and increases as 
the system becomes more transparent {R^{ec — 1) — ?> 0), while the skew-symmetric modes occur at a discrete set of 
values of h in the whole range h G {R, oo). 

The existence of the above classification can be established by studying the symmetry of the function e{x,z). 
Consider the operator Pa defined by PaEi^{x,z) = e~'^°''^''Ei^{x + a,—z). For a = or a = 5, the operator Pa 



commutes with the operator H of Eq.© and has eigenvalues ±1 since it is a projection, i.e., P^ = 1. By virtue of 
the commutativity, each bound state E^^ is an eigenfunction of Pa and, hence, 

E^{x + a, -z) = ±e'°'=-£;^(x, z) . 

Symmetric states may then be defined as those for which E^{x + a, — z) = e'^"'''^E^{x,z) while the skew-symmetric 
states satisfy E^{x + a, — z) — —e^°'^''Et^{x^ z). For a = this simply means that the first states are even in z while 
the second are odd. 

This description cannot be generalized to an arbitrary shift a since Pa no longer commutes with H if a G (0, §). 
Instead, it follows from the symmetry of £{x,z) that the operator QaE^{x, z) ~ E^{a — x,—z) commutes with H 
for all a and also Q^ = 1. However, the operator Qa is not symmetric in the subspace of functions in L2{S^ x R) 
satisfying condition ([3]) unless fc^ = or fc^; = ±7r. Indeed if it were, then any bound state Ei^ would have been an 
eigenfunction of Qa, and, therefore, Ei^{a — x, ~z) — ±Ei^{x, z). By replacing x by x + 1 and applying condition ([3]) 
one gets e^'^'-'E^^x, z) = e~^^''E^{x, z) so that fc^^ = or fc^; = ±7r. 

As noted above, the bound states are found by perturbation theory in the limit of thin cylinders. The technical 
details are given in Appendix [B] where it is shown that non-trivial solutions to Eq.([7]) can only exist if the fields 
Eu]{—he^) and Ei^{aei -\- he^) on the two cylinders positioned at (0, —h) and (a, h) respectively do not simultaneously 
vanish and satisfy the homogeneous system of equations 



(9) 



I <^oEUaei + he-i) + ^+EU-he-i) = 
\<^-E^{aei + he^) + '^„E^{-he^) = 

in which the coefhcients $o and $^ are given by the following expressions 

$o(fc,fc,)= Y. fr^-T-nA^rrVf (77^ + 2 in(27ri?) 

mf^ooV^^,™ 27ri(|m|-|-l)y 2t: \5o{k) 

°° i{±a{k^+2Tmi)+2hk^ „) 

<D±(a,/i,fc,fc,)= ^ — (10) 



with kz,m = \/k^ — {kx + 27rm,)2 being the z-component of the wave vector in the m diffraction channel and 

So{k)^(^) (£c-l) (11) 



In particular 5o{k) > as Sc > 1- It should also be kept in mind throughout the work that So{k) ^ 1 in the limit 
considered kR ^ 1. In terms of the fields on the cylinders in Eqs.®, the electric field strength everywhere off the 
scatterers is then given by, 

/ i((£c-a)(fc^+27rm) + |z-h|fe^,„) i(x(k^+2TTm) + \z+h\k^,m)\ 

E^r) = 2TTi6oik) E^aei + he^) ^ + E^i-hes) ^ (12) 

V ™ '^z.m '^z.m I 

The system of Eqs.Q admits non-trivial solutions if and only if its determinant 

A(a,/i,fc,fc,) = $2_^+^- (j3) 

vanishes at some point (a, h, k, k^) in the space of system parameters. This is the condition for bound states to exist, 
no matter if they are below or in the radiation continuum. 

In the following two subsections, roots of A(a, h, fc, kx) in k for fixed a,h, and kx are analyzed to find bound states 
below the radiation continuum as well as bound states in the continuum I. The analysis of the higher continua in 
the spectrum, while being similar to the case of the continuum I, is technically more involved. To avoid excessive 
technicalities before the discussion of applications of resonances with the vanishing width to a near field amplification, 
the analysis of higher continua is postponed to Section IIVI In each study, the spectral parameter k ranges over an 
open interval in which the functions $0 and $^ are analytic in k and diverge at the endpoints. Also, note that, since 
A(a, /i, k, kx) is even in fe^, the range of this parameter can be reduced from [— tt, tt] to [0, tt]. Throughout the rest of 
this section as well as in Section ITVl the symbols $*, $£, and $s denote, respectively, the following functions: 

cos(27ram) , $s(a, h, k, kx) — y. sin(27raTO) 

, Qz,m , Qz,m, 

(14) 



8 

where the superscript "cZ" in rrf^ means that the sums are taken over all m's that correspond to closed diffractio n 
channels, and qz,m is the imaginary part of k^.m when the m*^ channel is closed, i.e., qz^m = \/{kx + 2nmy — k^ if 
fc^ < {kx^-2TTmY. Recall that kz^m = \Jk'^ ~ (kx + l-nm)'^ . By construction, the functions P^ are always real- valued. 

A. Bound states below the radiation continuum 



Bound states below the radiation continuum are nontrivial solutions to the homogeneous Lippmann-Schwinger 
integral equation when all diffraction channels are closed, k^ G (0, vr] and < fc < fc^;. In this case the determinant of 
Ea.([T^ factorizes as 



A(a,/i,fc,fc:r) = -*+*", *± = $* T ^J^l + *? 



(15) 



hence bound states will exist for wavenumbers k at which the functions ^I'^ and ^ vanish. That such wavenumbers 
exist for ^+ follows from the limits 



lim *+(fc) 

fc-i-0+ 



-oo , 



lim *+(fc) 



and the Intermediate Value Theorem. In particular, for each triplet (a, /i, fc^;), there exists a wavenumber k^{a, h, kx) 
at which ^+(fc+) = 0. It is shown in Appendix [Dl that in the leading order of 5o{kx)^ 



kx 






(16) 



The function ^~ can also be shown to have roots when the distance 2h is sufficiently large or when kx is close to 
TT, this gives a second family of wavenumbers k~ (a, h, kx) at which bound states occur. These assertions result from 
the limits at and at kx of ^^. These limits are, 



lim ^ (fc) = oo , 



lim vj/ (fc) = -4/i cos (ttc) 



sin (tto) 
^7r(7r - k^ 



+ 



1 



2Ti5„{kx) 



0(1), 



In particular for fixed kx G (0,7r] and a j^ ^, h can be chosen sufficiently large for the last limit to be negative and 
therefore by the Intermediate Value Theorem a wavenumber k^ such that ^/^{kr) ~ exists on the interval (0, kx). 
Conversely, if h is fixed, then 



lim lim ^ 



-Ah 



2'k5o{ti] 



if a 7^0 
+ 0(1) ifa = 



(17) 



Thus if < a < i and kx is sufficiently close to tt then the equation '^^{k) — has a root k^{a, h, kx). When a — 0, 
the sign of the limit in Eg. dTTI) depends on the radius of the scattering cylinders, the dielectric constant Sc and the 
distance between the two arrays. If h is large enough or R^{£c — 1) is not too small, then the limit is negative and 
the existence of the wavenumber k~ follows. Because of a complicated dependence of the existence condition for the 
second family of wavenumbers on the physical parameters of the system, such a simple analytic approximation as for 
k'^ above is not possible for k^ . 

Finally, from Eqs.® and ([T2l) the bound state solutions E^ at the wavenumbers fc* are obtained. As any solution 
to a homogeneous equation, they can only be determined up to a normalization constant which is chosen to be the 
value of the electric field E'^ (—he^) on the cylinder at (0, —h). In terms of this value, the electric field on the cylinder 
at (a, h) is then 



E^iaei + tea) = ±e'^'^+'''''^^E^{-he3), 
and everywhere off the scattering cylinders it is: 



E^ 



2m5o{k)E^{-hez 



E 



gia;(fca;+2irm) 



arg 



^i\z+h\k^ 



E 






k^k± 



±e 



i{\z — h\kz^m — 27Ta7n-\-(f)) 



(18) 



k=k± 



B. Bound states in the radiation continuum I: One open diffraction channel 

When only the 0"^-order diffraction channel is open, k^ G [0,7r) and kx < k < 2tt ~ kx, the determinant of Eo. lfTS)) 
can be rewritten in the following form convenient for the analysis 

A(a, h, k, kx) = £i!!!i|^ + $2 _ VI/+VI/- + |1 U,+ sin^ihk,,) + *- cos^hk,)) (19) 

fc^ kz 

where '^^{a,h,k,kx) = ^^{k,kx)± ^EHi^Mi) _ ^^(a,h,k,kx)j and the fmictions $*,$c and $s are given in Eqs. (fT4) ). 
Thus, bound state exist if both the real and imaginary parts of (J19p vanish: 

sin^(2fefc,) ,^ 

fc? '■ (20) 

*+sin2(/ifc^) + *-cos2(ft.fc^) =0 

The first of these equations implies that ^Jz+v]/^ > 0. If this inequality were to be strict, then the second equation 
would not have held, and, therefore, \I/+\I>~ = 0. Thus, the first equation implies that $« = and sm{2hkz) = 0. 
In turn, the latter equation implies that either cos{hkz) = or sm{hkz) = 0, and, therefore, the system of Egs- ipO)) 
splits into two systems, namely, 

<^s{a,h,k,kx) = ((^s(a,h,k,kx)=0 

cos{hkz) = < sm{hk,) = (21) 

'i+{a,h,k,kx) ^0 [*"(a,/i, fc, fcx) = 

To solve the first equation of each system, the series for "I>s is rewritten as, 

$s = - 2J Cm sin(27ram), c„i = (22) 

_, Qz^ — m Qz.ni 



Recall that qz^m = \/{kx + 2nmY — fc^- In particular, if fca; = then Cm — 0, Vm — 1, 2, 3... and the equation $s = 
holds trivially. Similarly, if a = or a = i, then the equation holds trivially as sin(27ra7n) — 0,Vm. It turns out that 
these are the only possible roots of $s if ft, > ^ sa 0.055. This conclusion stems from the following factorization of 



$s(a, ft. A;, kx) — — sin(27ra) y 

m— 1 

together with the facts that for fc^, 7^ 0, 



Cm — 2Cm+l + 2 y ^ {Cm+2n — Cm+2n+l) 



sin (Tram) 
sin^(7ra) 



(23) 



Cm>0 and sup i ^^^^ I = e-4^" (24) 

m=l [^ Cm J 

These statements are established in Appendix [Cl When ft > ^, then Cm+i < ^Cmyrn and therefore each of the 
summands in the series of Eg. (1231) is nonnegative. Since the first term of the said series does not vanish, it follows 
that the series does not vanish. Consequently, if fc^, 7^ then <I>s = if and only if sin(27ra) = ,i.e., a = or a = |. 
It is remarkable that under the restriction ft > ^^ no solution to Eas. ((2T|) is lost. Indeed, since sm{2hkz) = at a 
bound state, it follows that 



for some positive integer n. In the ranges considered, kz < 27r and therefore h > j > j > ^^. Thus a necessary 
condition for the existence of bound states is that either kx = while a e [0, i] is arbitrary or a is either or i while 
kx (z [0,7r) is arbitrary. In set notation, 

(a, kx)eL^ ([0, i] X {0}) U ({0, 1} X [0, vr)) (25) 
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The set L is represented by the second level of Figl2][a) . 

Let us turn to solving the last two equations in each of the systems in Eas. ipT]) . For this purpose, the function ^^ 
is defined for each positive integer n by, 



'i'^(a,-——,k,kx) if n is odd 
2fez 



("'2fc~' '^^ it n IS even 



^nik,ka:,a) 



27r<5o(fc) ^^ \ 27r(|m| + 1 

where k G {k^, 2n — kx). Then the systems of Eqs. ([2T|) split into the countable set of systems 

m: 




2y/k^ - kl n = l,2,3... (27) 

_*„(A:,fc^,a) =0 

where the systems corresponding to odd n arise from the first of systems (j2ip and those corresponding to even n 
result from the second system. In each of the systems p7p . the second equation determines the wavenumbers at which 
bound states occur. In turn, by the first equation, these wavenumbers determine the distances h that allow for the 
bound states to exist. In the next paragraph it is shown that, for fixed {a^k^) £ L, each system can admit at most 
one solution. Hence the set of distances h allowing the existence of bound states is discrete. 

Let (a, kx) G i be fixed. It is shown in Appendix [Q that the function k i— >■ '^^ is monotone decreasing on its domain 
(kx,2TT — kx) and therefore admits atmost one root in the said domain. Moreover, the root only exists if the limits 
of fc I— >■ ^„ at kx and 27r — kx are of opposite sign. Specifically, the limit at kx must be positive while the limit at 
27r — kx must be negative. The first limit is, 

<x> if fcj, = 
lim *„ = <j 1^ y- (^^ ^_^_^^^_^+i('i+ln(27ri?)^ if fc, ^ 



k^kz \ 2'K5o{kx) ' ^V27r(|m| + 1) ^^Tr'^m^ + Aurnkx ) ' tt V2 



The requirement that this limit be positive when kx ^ Q puts a restriction on the values of R and Sc that allow for the 
existence of bound states. However, this is a too complicated condition to analyze. A weaker, but easier to analyze, 
condition is obtained by first rewriting the positivity condition as. 



•^ f 1 1 I \ I fl 



m— 1 



> >. .... . 7 + .. . .7 . . 7 -— +zU-M27ri?) (28) 



7ri?2(e^ - l)fc2 ^ y^/Air'^m? - ATrmkx ^/Aii'^m? + 4nmkx ti"™/ ti" V^ 



The rearrangement of the series is made to ensure that all the summands in the series are nonnegative; they vanish 
at kx = 0. Also, since the cylinders are thin, it may be assumed that R < 2n ~ 0.262 so that -i — ln(27ri?) > 0. 
Therefore all the summands in Eq. (pS)) are nonnegative and hence the left hand side must be larger than each of the 
summands on the right individually. To estimate the threshold value of Ry/sc — 1 below which bound states may 
exist, the first term {m = 1) is retained in the sum (j28|) . This term is then written in the form A:^(7r — kx)~^''^g{kx) 
to isolate its branch point at fc^, = tt and its root at fca, = which is of multiplicity 2. The minimization of g{kx) on 
[0,7r] produces an estimate: 



C^/¥^kx 



RV^^^T < ^-72 — - (29) 

where 



C^.t72 min ^ ^ l + VTTt+yr-t ^^_^^^ 

\ o<t<i ym (1 + vi-t^) (2 + vT+t + VT^) J 



to observe In particular when kx is close to tt, the quantity Ry/Sc — 1 must be small enough in order for bound states 
to exist at all. Figure [5Jb) shows the regions in which Eq. l^^ is valid. 
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As specified above, for bound states to exist it is required that the hmit of '^n at 27r — k^ be negative. This Umit is, 

— oo if n is odd and a ^ ^ or n is even and a ^ 



hm *^ 

fe-J-(27r-fc^)- 



^7r(7r - k^) 2Trdo{2Tr ^ k^) 



0(1) otherwise 



In particular, the limit is negative except possibly when n is odd and a = i or n is even and a = 0. Even in the 
latter cases however, the negativity condition may be ensured by taking k^ sufficiently close to tt so that k~^ is large 
or by choosing n sufficiently large. Thus, if the parameters R and £c of the scattering cylinders verify condition (|28p . 
bound states in the continuum I do exist. To be precise, given a positive integer n; then 

• y{a,kx) e ((0, |)x{0})u({0}x [0,7r) j there exists abound state at the wavenumber fc2„--i (a, fej^) € (fcx, 27r— fcj^) 
and at the distance 2h2n^i{a,kx) ~ j- " '"^ ^ between the two arrays of cylinders. When a — ^, then the 
wavenumber fc2n-i exists for sufficiently large n or for k^ sufficiently close to tt. 

• V(a, kx) G ( (0, |) X {0}) U ( {^} X (0, tt) j there exists a bound state at the wavenumber fc2n(a, k^) G {k^, 27r~-kx) 
and at the distance 2/i2„(a, fcj^) = , ^"^ between the two arrays of cylinders. When a — 0, then the 
wavenumber fc2„ exists for sufficiently large n or for kx sufficiently close to tt. 

In the limit of the thin cylinders considered, approximate values can be inferred for the wavenumbers fc„, n — 1,2,3.., 
by only keeping the leading terms in the equations ^n(fc) = 0. This is detailed in Appendix [Dl For instance if n is 
odd and a = or ri is even and a = |, then the wavenumber knikx) and the distance hn{kx) are approximated in the 
leading order of ^o by, 

MM^^,-..- ^-^^-^-' ..(..). ^4:L=(i + ^"-'^7 ^-' ) (30) 

Svr - kx 4^7r(7r - fc^;) \ tt ~ kx J 

Finally, from Eqs.® and ((T2|) the explicit form the electric field for the bound states {-Enj^i is obtained at the 
wavenumbers {fc„}5^i and the distances {2/i„}^j^ between the arrays of cylinders. The eigenfunction of bound states 
in the continuum can only be determined up to a multiplicative constant which is chosen to be the value of the electric 
field £'„(— /inGa) on the cylinder at (0, —hn). In terms of this value, the electric field on the cylinder at (a, /i„) is then, 

and everywhere off the scattering cyhnders it is: 

i?„(r) = 27rz5o(fc„)i?„(-/i„e3) ^ -^^ (^e^N+'^-lfe?.™ + (_i)«+ig,(|.-^„|fer,„-2.a™) j (3^^ 



where fc"^ = ^fc^ — (kx + 2'nmY . It can be verified easily that outside the scattering region, i.e., \z\ > ft.„, there is 
no contribution to the field E^ from the 0— order diffraction channel. With this channel being the only open channel, 
it follows that En decays exponentially in the asymptotic region |z| — >■ oo and therefore it is square integrable on 
S^ X R a,s required. Figure |3] shows examples of plots of the absolute values of the fields £"„. 

C. Application: Zero width resonances and near field amplification 



The study of Section Hill shows that, given two parallel arrays of subwavelength dielectric cylinders of radius R and 
dielectric constant Ec, there are points (a, kx) for which a bound state exists if the distance between the arrays attains 
a specific value 2/if,; this value also determines the wavenumber fc(, of the bound state. Throughout the following 
discussion, the pair {a,kx) is fixed. 

Consider the scattering problem for the double array when h — h^. If the incident wave has the wavenumber fch, 
then the solution to Eq.Q is not unique as any solution of the homogeneous equation can always be added. The 
latter are the bound states. Of course, this ambiguity is related to the fact that the incident radiation cannot excite 
the bound state (which is a wave-guiding mode propagating along the array), and, hence, an addition condition must 
be imposed that the bound state is not present in the system at the very beginning if one wishes to have a unique 
solution. On the other hand, the presence of a bound state has no effect whatsoever on the scattering amplitudes 
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FIG. 3. (Color online) The modes E„ as defined in Ea. (|3ip for a — 0, R = 0.1, Ec ~ 1.5, and kx — 0. The panels show En as 
a function of x (vertical axis) and z (horizontal axis). The color shows the absolute value of E„ as indicated on the left inset. 
The positions of cylinders are indicated by a solid black curve. The values of /i„ are shown below each panel. Top left: Ei, 
the symmetric mode for the smallest h = hi. Top right: E2, the skew-symmetric mode at h — h2. Bottom left: E3, the second 
symmetric mode at h = ha. Bottom right: E4,, the second skew-symmetric mode at h = hi. 



as they are defined in the far field zone (z — >■ ±cxd) to which the bound state gives no contribution anyvifay. This 
degeneracy disappears as soon as {h,k) ^ {hb,ki,). In the latter case, the solution E^ to Maxwell's equations is 
uniquely determined by Eq.®. Such irregularity suggests that, as a function of h and k, the field E^ is not analytic 
in the vicinity of the points {hb, kb). This is indeed the case. It is shown in what follows shortly that the values of the 
reflected flux as well as those of the fields inside the cylinders near the points (/ib, fcf,) depend on the path along which 
these points are approached in the h, fc— plane. From a mathematical point view, this lack of analyticity is explained 
by the presence of simple poles at the wavenumbers fcf, in the field E^ when considered as a function of k for fixed 
h = hi,. The objective here is to exploit the existence of these poles to show that the evanescent field in the scattering 
problem is amplified as compared to the amplitude of the incident field when (/i, k) is close to (/if,, fct) in some regions 
of the array, in particular, on the cylinders. The effect can therefore be used to amplify optical non-linear effects in 
the structure in a controllable way. This is illustrated with an example of one open channel for an array without the 
shift, i.e., a = 0. 

Suppose that a plane wave of wavenumber k € {kx, 2Tr — kx) and unit amplitude impinges the double array. In this 
case, the specular reflection coefficient which is the ratio of the reflected flux to the incident flux at the spatial inflnity 
is, 

7^ = |i?„p 

where Rg is the reflection coefficient of the only open diffraction channel, namely, the order channel as given in 
Eq. (|5b[) . If {h, fc) is not one of the points (hi,, fc&), then the refiection coefficient Rq and the fields inside the cylinders 
are, 



i?o — — 



cos^ (hkz) 



cos-^ {hkz 



E^{±he3 



ikz 



\ikz-^+ 

COs(/lfcz) 



sin {hkz 



sin {hkz) 



±i- 



\ikz-^- 
sin(ft,fcz 



27r(5o(fc) \cos'^{hkz) + \ikz'^+ sir? {hk z) + kikz"^ 



(32a) 
(32b) 



for the functions ^* of Eq.lfTOl) (See Appendix [B] for details on the derivations of the above expressions). The 
denominators in both expressions are factors of the determinant A(0, /i, fca;, fc)in Eq. ([T9| (Here a = 0). In particular, 
the points {h^, fc;,) are roots of the denominators to the specular coefficient and the fields inside the cylinders. 

For the illustration purpose, the behavior of the specular coefficient and the fields inside the cylinders are studied as 
(/i, fc) approaches a critical point (/i^, fc^) = (ft,2n-i, ^2^-1) for some positive integer n in the ft., fc— plane. As described 
in Subsection nil B[ at the point (ft2n-ij fc2n-i) the following system holds, 

{cos,{hkz) = 
*+(/i,fc)=0 

In the rest, the curves cos{hkz) = and '^^{h, fc) = will be denoted by Cc and C"*" respectively and their intersection 
points, i.e., the points (/i2n-i, fc2n-i), will be denoted by P2n~i- 
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The first observation is that as {h, k) approaches P2„_i, then cos{hkz) —>■ and ^~^{h, fc) — >■ independently as the 
curves Cc and C"*" intersect at a nonzero angle at P2n-i- This may be established through the linearizations of the 
functions {h, k) h-)- cos{hkz) and {h, k) i-^- kz'^^{h, k) at (/i2n-i, fc2n-i)- If Aft, = h — /i2n-i and Afc — k— fc2n-i, then 
in the vicinity of P2»i-i, 

cos(/ifc.) « ^(Aft, Afc) = (-1)" j^ ^2"-ife2"-i ^fc ^ ^^ 2„-l A/i 

V fcz.2n-l 

11/ 

-fc^^+(ft, fc) « 7j{Ah, Afc) = -fc^,2«-l (^afc*+(ft2n-l, fc2n-l)Afc + 9,,*+ (ftan-l , fc2n-l ) A/l 



where fc2,2n-i = \/^2n-i ^ ^i- The functions ^ and r/ are then linearly independent if, 

'-p ^^-— 9/i*^(ft2n-l, fc2n-l) - kz,2„^idk'^^ {h2n-l, k2n-l) 7^ 

fcz,2n-l 

That this condition indeed holds can be proved by examining the function ^2n-i(fc) = ^^ ( '^fc ; ^ ) introduced 
in Eg. ipS)) . In Appendix [Cj it is shown that dk'^2n~i{k) < for all fc e (fc^^, 27r — k^)- Consequently, 

5fc^2n-l(fe„-l) - -%^i^^^^9/.*+(/i2„-l,fc2„-l) +afc*+(/l2„-l,fc2„-l) < 

This establishes the linear independence of ^ and ry. Thus, as (A/i, Afc) — > (0,0), there should be ^ — ?► and 77 --> 
independently. In the vicinity of the critical point P2n-i the principal parts of R„ and Ei^{±he3) are then, 

R(hk)- ^ I eiAh,Ak) 

"°^"' ^ -^ 1 + ifc,,2„-iVl/-(/i2„-i, fc2„-i) ^^(Aft, Afc) + ^r^{Ah, Afc) ^^^''^ 

ikz,2n-i f . (-1)"+! ^(Aft,Afc) \ 

"^^ '"'^ "" 2^5o(fc2„-i) l^*l + zfc.,2n-i*-(/^2„-i,fc2„-i) ^e'(Aft,Afc)+*,7(Aft,Afc)j ^ ^ 

The first summands in each of these equations are constant and obey the estimate, 

1 So{k2n-l) 



l + ikz,2n-l^ {h2n-l,k2n-l) 



^z,2n-l 



The second summands in Eqs. (|33l) account for the lack of analyticity of the specular coefhcient TZ and the fields 
E^{±he3) in the vicinity of the critical point P2n~i- In particular, since ry = along the tangent line to C"*" at P2„_i, 
it follows that along this tangent line and hence along the curve C+, 

E^{±he3) « , ' , TTTX ' Ah = h- h2n-l -> 



2-<5„(fcL_,)(l-^i2f^i^||l)Aft' 



Thus the electric field inside the cylinders diverges at the points P2n-i- 

For the specular coefficient, the significance of the points P2n-i is that they are positions of resonances with the 
vanishing width, which, in turn, demonstrates that the bound state in the radiation continuum are interpreted as 
resonances with the vanishing width in the formal scattering theory. Indeed, if ft 7^ ft'2n-i is fixed and fci.(ft) is a 
wavenumber such that ^"^(ft, kr{h)) — 0, then fc^ is a resonant wavenumber for the specular coefficient 7?.. For fc near 
fcr, the Breit-Wigner theory asserts that Rg will have the form, 

R„ 



f\i tZj- ~p tJ. 



where 2r is the resonance width of the Lorentzian profile of 7?. = |PoP- The half- width F may be found by expanding 
the function fc 1— > ^"'"(ft, fc) in a Taylor series at the resonant wavenumber fc^, it is, 

2cos^(ftfc2 



kzdk^+ 



k=kr. 



so that at the points P2n-i this width vanishes. Figure U shows plots of the specular coefficients and the electric field 
along the curve C+. 
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FIG. 4. (Color online) The specular coefficient and the electric field on the cylinders near a bound state in the continuum for 

a = Q, R = 0.1, and Sc ~ 1.5. 

Panel (a): Shows the specular reflection coefficient as a function of h and the wavenumber k. It is plotted for kx — j near 

the threshold fc_i — ^. The specular coefficient is very close to its maximum along the curves <& (ft, fc) = that determine 

the resonance positions. The lower (upper) bright region roughly corresponds to the curve '^^{h,k) — {'^~{h,k) = 0). The 

bound states correspond to the points separating two consecutive bright strips. At the bound states the specular coefficient 

lacks analyticity and its value depends strongly on how the bound state is approached in the {h, fc)-plane. 

Panel (b): The specular coefficient TZ{k{h), h) (dashed blue curve) and the absolute value of the electric field Eui{he3){k{h) , h) 

(solid red curve) on the cylinder at {x,z) — (0, /i) of the double array where k — k{h) is implicitly defined by ^^{h,k) — 0. 

Along this curve, the electric field diverges near the bound states. 

Panel (c): Same legend as for Panel (b) in the case of two open channels and kx — it. In this case too, the fields on the cylinders 

diverge in the vicinity of the bound states. 



IV. BOUND STATES IN THE RADIATION CONTINUUM N, N> 2 



When more than one diffraction channel are open, the bound states may still be shown to exist. However they 
become rarer as the number of diffraction channels increases. The physical reason for that is simple. As suggested 
in the introduction, a bound state is formed due to a destructive interference of the decay radiation of two quasi- 
stationary electromagnetic modes localized in the vicinity of each array. If more than one decay channel are open 
for these modes, then the destructive interference must occur in all the decay channels in order for a bound state to 
form, which puts more restrictions on the system parameters. Indeed, consider, for instance, possible choices of the 
parameters a and k^ that allow for the existence of bound states. It was shown that bound states do exist below the 
radiation continuum for all pairs (a, kx) ,i.e., no restrictions at all. When one diffraction channel is open, then bound 
states can form when the pairs {a,kx) lie on the set L defined in Ea. ([25t . When two diffraction channels are open, 
then the bound states are shown to only occur if the shift a is or | and for a specific dense set of values of kx ■ As 
one goes on to higher levels of the spectral cylinder in Figlljd), the values of kx at which bound states may exist 
become more sparse and are determined by solutions of a system of diophantine equations. The conditions under 
which bound states exist are formulated first for the studied system in the continuum N > 2. 

The above assertions follow from the observations that, if the diffraction channels 0,-1,1,... are open, then at a 
bound state, the parameters kx , k , h ^ a , R , and Ec satisfy the relations. 



'A(a,/i,fc„fc)=0 
sm{2hk^) = 
8m{2hkz,-i) — 
sin(2ftfc34) = 



(34) 



where A(a, /i, fc^,, fc) is the determinant of Eq. ([l3|) . The additional equations are a result of the square integrability of 
the bound states on S^ x M. Indeed, if a solution Ei^ to Eq.([8]) is to be square integrable on S*^ x R then, the function 



z i-^ 



\E^{x,z)fdx 
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is integrable in z over R. Now, 



'y— 



I y^ 1 

^vi o-p Z , 7M 



ihkz 



-ihkz 






z -^ — oo 



where the superscript "op" in to°^ indicates that the summations are to be carried over all m's that correspond to 
open diffraction channels. In particular, for square integrability to hold, each summand in the equations above must 
be zero. It follows that for each open channel TOq the following system holds. 






(35) 



For nontrivial solutions, the determinant of the above system must be zero. Thus, siii{2hkz^mo) = for each open 
channel mg. Moreover, by considering the ratio of the field E^{aei + he^) to the field E^i^—he^) for all the open 
channels it is deduced from system ([55]) that. 



cos{2hk,) = cos{2hk,^^i)e-^''"' = cos{2hk,^i)e'^ 



(36) 



In particular, since for each open channel nio we have sm{2hk2,^mo) = 0, then cos{2hkz^mo) — il- Thus e^'^*" — ±1, 
hence a = or a = |. 

The first two equations of system (I34p determine k and h as functions of k^ while the last equations determine the 
values of fc^;. Thus as pointed out above, the values of k^ at which bound states may exist become more sparse as the 
number of open diffraction channels increases. 

In the case of two open diffraction channels, it is remarkable that these values of kx are dense in [0, tt]. A proof of 
this statement is given in Section IIVAI It is shown there that the values of kx allowing for the existence of bound 
states occur in a double sequence fc"'' where n, I are positive integers. In the leading order of R'^{ec — l), the elements 
of the subsequence fc^""*"^'' are shown to have the form : 



„2n+l,l 



2r2-l 



7r5(r2 - l)(4r2 - 1)'^ 
4(2r2 - 1)5 



RHsc - lY 



I 



2n+ 1 



(37) 



The elements of the subsequence A:^"'' are harder to derive due to a more intricate dependence on the system param- 
eters. For the subsequence fc^""'""'^' , the corresponding wavenumbers and distances between the arrays at which the 
bound states occur are proved to be obtained by substituting A;^"+^'' into the following expressions: 



k2n+i{kx) ~ 2n + kx - 



8Tr^S^^i2n + kx) 
2'K + kx 






^,52(27r + fc,) 



(38) 



As in the case of bound states in the continuum I, bound states in the continuum II are shown to only occur under 
the following restriction on the the radius and dielectric constant of the scattering cylinders : 



R^c 



1 < 



\/tt - kx' 



C« 2.016 



(39) 



In particular, when k^ is close to 0, the quantity R\/£c — 1 must be small enough in order for bound states to exist 
at all. Figure [2][c) shows the regions in which Ea. (l39)) is valid. 

The near field amplification observed in the case of one open diffraction channel persists when two channels are 
open. This can be established via an analysis similar to that of Section IIIICI Figure HJ^c) gives an example of such 
an amplification. 

When three or more diffraction channels are open; the equations sin(2ft,fcz_m„) — for each open channel TOq, 
determine the parameters kx-, k and h. In fact, if no, ni, n2, ... are positive integers such that 2hkz = noTr, 2hkz,-i — 
niTT, 2hkz^i — n2TT... then 2nQ y^ nf + n^ and, 



kx — 



2nl 



n\ 



h^ 



4^/2 



2^0 — n\ — 712, k = 



\/{n\ + n\- Anly + An\n\ 



2nl 



nl 



(40) 
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When only three channels are open, the integers ng, rii and n2 are only required to be in the order uq > ni > n2. 
When four channels are open, the additional equation sm{2hkz^2) = in system ([M| requires that the aforementioned 
integers satisfy the system, 

J3nl +nl =3nl + nl 
[no > ni > n2 > n^ 

As more diffraction channels become available, there are more and more constraints on the integers Ui, i — 0,1, 2... 
Provided integers satisfying those constraints can be found, bound states will then be formed for double arrays for 
which the radius R and the dielectric constant Ec satisfy A{a,h,kx,k) — with h, k and kx given by Eqs. (|40l) and 
a e {0, |}. From Eqs. (ITU|) and (IT^ . it follows that R and Ec must be on curves, 

o 

ln(27ri?) = C(no,ni,...) 



where k is given in Eq.(|40]) and C is some constant that depends on the integers n^, i = 0, 1, 2, ... 

A. Bound states in the radiation continuum II: Two open diffraction channels 

Suppose that both the O"^ and —1''' diffraction channels are open ,i.e., k^ € (0, tt] and 2tt ~ k^ < k < 2Tr + kx- 
Conditions (IMl) then translate to the following system of equations: 

2(1 — cos(27ra) cos{2hkz) cos{2hkz^-i)) t + t - 

^ l' 1 — cos{2hkz) 1 — cos(27ra) cos(2/i/i;2^_i)\ Z' 1 + cos{2hkz) 1 + cos(27ra) cos(2/i/i;2__i)\ 

kz fc^._i J \ kz kz.-i J 

(41) 
where ^^ = <i>* =F <i>c for the functions $* and $c of Eqs. p3) . Also, recall that o is necessarily or i as derived from 
Eqs.dSH). 

The first equation of the system implies that ^1'+^^ > 0. If this inequality were to be strict, the second equation 
would not have held, and, therefore, vl/+v[/~ — Q. Thus either ^+ = or 'f^ = 0. Note that the functions ^1'+ and 
^~ cannot vanish simultaneously. This can be verified by observing that, 

$+2 + ^-2 = 2{<i>l + ^l) 

and therefore, if the functions ^+ and ^^ were to vanish simultaneously; it would follow that $c = 0. But, 

22 > if a = 

™#o.-i ^^'^ 

Y^ (-1)™ < ff a = i and fc^ 7^ vr 

, Qz,m 

if a = 2 and kx = tt 

Thus (a, kx) = (|, vr). But then kz = fcz.-i and the first equation in system (I4ip reads, 

l + cos^{2hkz) = 

This is impossible and therefore at a bound state ^+ and ^P^ do not vanish simultaneously. In particular, there are 
no bound states corresponding to the pair {a,kx) — (|,7r). This is the reason this point was removed from the third 
level of Fig. Ela). 

For the remainder of the discussion it is assumed that {a,kx) ^ (f,'?''). System (|4T|) then splits into two systems, 
namely, 

cos{2hkz) = —1 \ coa{2hkz) — 1 

cos(27ra)cos(2/ifcz _i) = -1 < cos(27ra) cos(2/ifcz _i) = 1 (42) 

*+ = U'" = 



$n= <^ 
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Thus to each value of a corresponds a pair of systems whose solutions, if any, give rise to bound states in the continuum. 
For a — Q, these systems are, 

{cos(2/ifc^) = -1 {cos{2hk^) = l 

cos(2/ifc^^_i) = -1 (B): \ cos{2hk^^^i) = 1 (43a) 

while for a = ^ they are, 

{cos{2hk^) = -1 {cos{2hkz) = 1 

cos(2/ifc2,_i) = 1 (D): \ cos(2/ifc^,_i) = -1 (43b) 

The existence of solutions to the last two equations in each system is proved first. Then the first equation of each 
system is added to show the existence of bound states. 
For each positive integer n, the function ^„ is defined by, 

TITT 1 

^^(a, — , fc, kx) if n is odd and a = or n is even and a = — 

^n{k,kx,a)^ { ^ ^^^ ^ 

^^(a, — , fc, kx) if n is even and a = or n is odd and 0=7: r..\ 

2/^2.-1 2 v^^i 

1 V^ / 1 l_(_l)"+2a(m+l)g-n^9.,„fc-Li\ W3 



y ^^ — -.-—^^^ — ^ — - — ——\+-[t+H2nR) 



27r5„(fc) f-^ ^ \ 27r(|m 



where k e {2ti — k^, 2it + k^)- Then the systems formed by the last two equations of each of systems (P5|) split into 
the countable set of systems 



2 VP _ (2^ _ ;,^)2 n=l,2,3... (45) 

\E'„(fc,fc^,a) =0 

It is shown in Appendix [C] that for each positive integer n, the function k i— > \l/„ is monotone decreasing on its 
domain (27r — kx^2'K + kx)- It follows that if system (j45|) has a solution, then this solution is unique. Moreover, such 
a solution will only exist if and only if the limit of fc i— > ^„ at 27r — kx is positive while the limit of ^„ at 27r + kx is 
negative. The first limit is, 

lim *„ ^ j^ ^ + V f —--^, ^ ) + i f ^ + ln(27ri?)^ (46) 

fc^(2.-fc.)+ " 2n5,{2n-kx) ^^_^\^27r(|m| + l) ^(27rm + fc,)2 - (27r - fc,)^^ tt V4 ^ 7 

As in the case of bound states in the continuum I, the requirement for this limit to be positive puts a restriction on 
the values of R and Ec that allow for the existence of bound states in the continuum II. An easily analyzable condition 
on these parameters is obtained by following the same procedure as in Section IlIIBl First, the positivity condition is 
rewritten as, 

2 ^ ^ / 1 ^ 1 1 \ 

TiR^iec - \){2tt - kxY ^ ^^ \ ^{2i:m + kxf-{2Ti-kxf ^ ^{2^{m + l) - kxf - {2tt - kxf TT^m{m+\)) 

+ i(,s-^-ln(27ri?)) (47) 



where 



f: ' 'f ' ' ' )V"- 



The rearrangement of the series is made to ensure that all the summands in the series of Ea. (|47)) are nonnegative; 



they vanish at kx = tt. Also since the cylinders are thin, it may be assumed that R < j-e'^ * w 0.150 so that 
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s — I — ln(27ri?) > 0. Thus all summands in Ea. (|47|) are nonnegative and hence the left hand side must be larger 
than each individual summand on the right hand side. To estimate the threshold value of Ry/Sc — 1 below which 
bound states may exist, the first term (m = 1) is retained in the sum (j47|) . This term is then written in the form 
(tt — kx)\/k^g{kx) to isolate its branch point at fcx = and its simple root at kx — tt. The minimization of g{kx) on 
[0,7r] produces the estimate ([M)) where C is exactly, 



C = 247r" 



(2-t)2 fVi^ 



Vt 



o<t<i V3~I \l + Vi V^+V3~t 



2.016 



As mentioned above, in addition to the requirement that the limit at 2t: — kx of fc i— > ^„ be positive, one must also 
require that the limit at 27r + kx be negative in order for system (j45p to have a solution. The latter limit is. 



lim *,, 



^/27rk'x 



0(1) 



^37r(7r - kx) 



if n is odd 

if n is even and a = 



+ 0(1) if n is even and a = ^ 



In particular, the limit is negative if n is odd. If n is even and a = 0, the negativity condition may be ensured by 
choosing kx sufficiently close to zero or by choosing n sufficiently large. If n is even and a = i, then kx must be very 
close to TT for the limit to be negative. For parameters R and Ec satisfying condition (|47)) . the conditions of existence 
of the solutions (kn{kx),hn{kx)) to system (|45p are summarized in Table ID In the leading order of 6o{27t + kx), 
approximate values of k2n+i and /i2n+i are given by Eq. psp fSee Appendix ID] for derivation). 

The approximate values for the wavenumbers k2n{kx) are more difficult to find due to the dependence of their 
existence on the physical parameters of the system. 





a = 


a-h 


n odd 


(k„, h„) exists yk^ G (0, tt] 


{kn,hn) exists Vfca; G (0, tt) 


n even 


(k„, h„) exists for n large or k^ small 


(fc„, h„) exists only for kx very close to tv 



TABLE I. Existence of solutions to systems (|45p . In particular, systems (A) and (D) in (|43p always have solutions whereas (B) 
and (C) might not. 



Since the conditions of existence of solutions to system (j45p are now established, the existence of solutions to 
systems (|^5)) can be investigated. So far only the last two equations in each of the latter systems have been used to 
determine the values kn{kx) and hn{kx) for each kx G (0,7r] that are susceptible to permit the existence of bound 
states. It follows that the first equations in each of systems (|33]) determine the values of kx that allow for the existence 
of bound states. In the coming paragraphs the set of these values is shown to be discrete and dense in [0, tt]. 

Let n be a positive integer for which kn{kx) exists for all kx G (0,7r] {kx G (0,7r) if a = i). Consider the function 
ifn defined by. 



ifinikx) = 2h„{kx)\/kl{kx) - kl = mi^ 



kn\kx) 



klikx) - {2n 



kx G (0,7r) 



where the value kx = t^ is purposely left out and will be discussed later. 

Collectively, the first equations of systems (j43|) are cos {ipn{kx)) = ±1 and therefore, they have solutions if the range 
of (fin can be shown to contain even and odd integer multiples of tt. That this is indeed the case follows from the 
continuity of ipn and its limits at and tt. These are. 



lim ipnikx) 



cx) and lim ipn{kx 



- niT 



so that the range of (pn contains the interval (n7r,oo). In particular, for each positive integer I > n, there exists a 
point kx''' G (0,7r) such that ipn{kx'^) = ^tt and therefore cos ((/5„(fc"'')) = (—1)'- This establishes the existence of 
bound states in the continuum II. As claimed, they exist for a discrete set of kx values in (0,7r), namely, the points 
A:""'. To each of these points corresponds a specific wavenumber kn^i = kn{kx'^) and a specific distance hn,i = hnikx'') 
at which a bound state in the continuum exists. Note that the points fc"'' depend on R, Sc and the shift a as is 
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illustrated for instance for the points /c^"+^'' in Ea. (|37|) (See Appendix ID] for details). Note also that by definition, 
the values fc2n+i,2;+i ^^^.^ solutions to system (A) of Eqs.(|33]) while the values fc^"+^'^' are solutions to system (B) of 
the same set of systems and that the two sets of points do not overlap. 

As far as the case fc^^ = tt is concerned, it was established in the beginning of this section that there can be no bound 
states if a = |. However, they do exist if a = 0. This is because kz — kz^-i a,t k^ = tt so that the first equations 
in systems (|43ap are superfluous and hence the existence of solutions to systems ([^S)) alone suffices to guarantee the 
existence of bound states in the continuum II (Note that systems (j43bp become inconsistent as expected). Thus if 
a = 0, the list of points {fc"'', / > n} is to be completed by adding to it the point kx = tt. Keeping with the notation, 
this point is A:"'" for each n since (fn can be extended to tt by defining (^„(7r) = mr so that cos((p„(7r)) = (—1)". Thus 
the point k^ — tt is of infinite multiplicity in the list {fc"'', I > n} and hence is associated with an infinite set of bound 
states. Indeed for each positive integer n, sufficiently large if it is even, there exists a bound state at the wavenumber 
kn,n — kn{'K) and at the distance 2ft,„,„ — 2/i„(7r) between the two arrays of cylinders. All the other points A:"'' change 
with the physical parameters of the cylinders and the indices n, I and therefore none of them is certain to be repeated. 

It is remarkable that the set of points fc"'' is dense in [0,7r] as is suggested by Eq. (|37| . To demonstrate this fact, 
consider an arbitrary interval I^^ = [a, (3) C (0, tt). It will be shown shortly that. 



lim 



i'P2n+l{a) - f27i+l{l3)] = oo (48) 



Therefore, however small the interval Iq.^ may be, the interval '^2n+i(Ia^) contains positive integer multiples of tt for 
sufficiently large n. If In is such a multiple, then fc^"+^'' e Iq^. Thus the points {fc^"^^''} are dense in (0, tt]. 

Before establishing the limit (US]), recall that fc„ = kn(kx) designates the solution to the equation '^n{kn{kx), kx,a) = 
0. The first task is to show that the sequence {k2n+i{kx)}^=i converges for each fixed kx- To this end, let ^oo be the 
function defined by, 

«'oo(A;,fc,) = -^--+ Y. ( o f\ V n '— V-fl+MS^^)) , k e {2TT - kx,2TT + kx) (49) 

For fixed (a, kx), the function k i— >■ ^oo is the pointwise limit of the sequence of functions k i— > ^„; it is also continuous 
and monotone decreasing on its domain (27r — /c^., 2TT + kx) as is shown in Appendix [Cl As fc — ^ {2tt + kx)~ , ^oo — ^ — oo 
while the limit of ^oo at {2tt — kx)~^ is exactly limit (^51) . In particular, the latter limit is positive as required by 
condition (|Tf)) . It follows that for each kx & (0, tt], there exists a unique point koo{kx) & {2tt — kx, 2tt + kx) such that 

^ cx>\l^cx>\l^x), kx) = U. 



Now for fixed kx. 



*oo(fc2«+l) = *oo(fc2«+l) - *2n+l(fc2«+l) 



OO 



s - 

m— 1 \ 



-(2n+l)7rq,^^fc-i_^ -(2n+l)7rq,^_^_iA;-i_^ ' 



oo 



Qz,m Qz. — m—l 



if a = 



^ ni—l 



^— (-1)"+! ifa=i 

Qz,7n Qz. — 711—1 J 

(50) 
where it is understood that in the two series on the right, the terms qz^^n and A:^ _i are to be evaluated at k — k2n+i- 
The sum of the first series is obviously positive. The sum of the second series is also nonnegative as it is the sum of 
an alternating series whose terms decrease in absolute value. Thus ^oo(A;2n+i) > 0. 
As fc iH- ^co is a decreasing function it follows that. 



'K2n+l 



(kx) < kooikx) < 2tt + kx , Vn = l,2,3... and Vfc^ S (0, tt] (51) 



In particular, Qz^i = \/{2tt + kx)"^ — k\^^^ > \J{2tt + kx)'^ — fc^ > and therefore qz^\ does not converge to as 
n —7- oo. Hence, 

/ -(2n+l)7rg,,ifc-Li\ 
*oo(A:2„+l) = O > = *oo(fcoo) 

Since the function k i— >■ ^I'co is continuous and bijective for each fixed kx € (0, tt], it follows that k2n+i{kx) — ^ kooikx) as 
n — )• oo. Note that without condition (ISTI) . qz^i could converge to causing the sequence {\E'oo(fc2n+i)}^i to diverge 
as indicated by Eas. lfSO)) . This is the reason the subsequence {fc2n}^i had to be excluded. For this subsequence, it 
can be shown through an analysis similar to the above that k2n{kx) ^ koo{kx) , Vn = 1, 2, 3... and therefore more work 
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would be needed to show that no subsequence of the sequence {^2™}^! converges to 27r + k^ 
From the convergence of the sequence {k2n+i{kx)}'^=i it is deduced that, 



1 / \ k^ (k ) — k^ 

In Appendix [C] the function (/?oo is shown to be strictly decreasing so that (/5oo(a) — '■PooiP) > 0. This establishes limit 
(H51) and thereby the density of the points {fc^""*"^''} in the interval (0,7r]. Moreover, if |A„(a,/?)| is the cardinality of 
the set An{a,P) = {I : fc"'' € (a,/3)} then, 



The last limit follows from the fact that as B?{ec — 1) — > then koo{kx) ^-271 + k^- 

Lastly, from Eqs.© and ((T2|) . the analytic expressions of the bound states £"„_; can be obtained for positive integers 
n and I such that I >n ii a — and I > n \i a = ^. As mentioned in prior sections, each of these states can only be 
determined up to a multiplicative constant which is chosen to be the value Enj^—hnje^) of the electric field on the 
cylinder at (0, —hn^i). In terms of this value, the electric field on the cylinder at (a, hn^i) is then. 



and everywhere off the scattering cylinders it is 

oi£c(fc"''+27rm) 



i?„,,(r) = 27ri<5o(A:„,0SrM(-/in,ie3) J^ -, [e\^+^r.AK:r. + (_l)n+2a(m+l)+lg.|.-/.„,,|fe.,™^ (53) 



^z^m 



where A:^;^ = ^^^^^ - (fcj'' + 27rm)2. 



Appendix 



Appendix A: The Lippmann-Schwinger integral equation 

In this section of the appendix, it is proved that the solution to the Lippmann-Schwinger integral equation in Eq.([6]) 
solves Eq.(l2]). This will be done by considering the functions involved as distributions acting on smooth functions of 
compact support. 

To do so, only locally integrable solutions to Eq.Q are sought. In this setting, all the functions involved in Eq.®, 
namely; e*"^"", E^j, eE^, {e — 1)E^, G and ((e — l)Ei^) * G, are locally integrable. 

Let then (^ be a smooth function of compact support and E^^ be a locally integrable solution to Eq. (JH) . Then 

{^E^ + eeE^,ip) = e{eE^,^) - e{E,,ip) + ^(((e - 1)^;^) * G, A(p) (Al) 

47r 

where Ei(v) = e'*'''. Since none of the distributions (e — \)Ei^ and G is compactly supported, the convolution used 
here is to be understood in the sense of the usual convolution of functions. Therefore, the last integral of Ea. (|Al[) 
may be interpreted as 



(((e - l)E^) * G, A^) = j (e(ro) - l)E^{ro){G,,,^^)dv, (A2) 

where Gro(r) = G(r|ro) and satisfies the distributional Helmholtz equation 

AGro + fc'Gro = -47r<5ro (A3) 

Therefore 

(Gro,A^) = (AGro,¥') = -A:2(Gro,^) -4^(5ro,^) 
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It follows that 

(((£ - 1)^^) * G, A^) = -k^ide - l)E^) *G,^)- ATT{{e - l)E^,ip) (A4) 

Thus (A_E^ + k'^eE^.ip) — for all test functions (p. 

Note that the only difference between the proof that Eq.® solves Eq.Q in the case of the finite array and that 
of an infinite array is the way Ea. (|A4p is derived from Ea. (|A2p . Indeed, if the array of cylinders is finite then the 
convolution of Eq. (jA2|) is in the distributional sense as (e — l)Ei^ would be a compactly supported distribution. 
Therefore one could establish Ea. (|A4[) immediately from Eq. (jA3|) and the identity of distributional convolution 

(((£ - \)E^) * G, A^) = (((£ - l)E^) * AG, ^> (A5) 

This identity is not immediate in the case of the infinite array because, as mentioned above, none of the convoluted 
functions has compact support. In fact, the function ((e — l)£'c^) * G is a conditionally convergent series so that 
integration by parts cannot be applied to establish Eq. (|A5p . 



Appendix B: Solution of the Lippmann-Schwinger integral equation: zero radius approximation 

Here an approximate solution to the Lippmann-Schwinger integral e quat ion is established in the small radius 
approximation. In the usual theory of scattering from small particles |19l - [23 ]. the problem is solved with high 
accuracy by assuming that far from the scattering region the solution is a linear superposition of the waves scattered 
by each individual particle. Due to the cylindrical geometry of the dielectric scatterers, it turns out that similar 
approximations may be made to solve the scattering problem on the double array but with solutions which are valid 
everywhere off the scatterers even in the region between the two grating structures. The validity of solution in this 
extended region can be used to find the so-called hot spots (see Fig. [3]) where the magnitude of the electromagnetic 
fields peaks. 

The structure considered is shown in Fig. [Ija). The cylinders in the structure are labeled as C„i^n where n is either 
1 or —1 depending on whether the cylinder is on the right or left array. The integer m refers to the x-coordinate of 
the cylinder's axis. In particular, for the right array cylinders, the axes are positioned at r^ i = (m + a)ei -t- ^183 and 
those of the left array are at positions r„j^_i = mei — he^^. 

The Lippmann-Schwinger integral equation may then be written as a sum over all cylinders as 

E^{v) = e^'^ + ^k\e,-l) ^ j E^^^^-^H„{k\r - r„|)dro (Bl) 

Far from the scatterers, each of the integrals is well approximated through the mean value theorem by 

/ EUro)G{r\ro)dro « in^ R^ e'"'""^ H„{k\r - r^.n\)EUro,n) 

so that the far field may be expressed in terms of the fields on the cylinders Go,±i as 

00 _. 

EUr)=e'''-'^ + nTS,^{k)J2EUro.n) ^ e'""=^i7„(fc|r - r„,„|) for 6„{k) = -k^R^e, - 1) (B2) 

n— ±1 rn— — 00 

The claim is that this approximation remains valid in the near region too. This may be established by means of 
the Bessel function expansions of the field inside the cylinders and the Hankel function Hg. To this end, let Tq be a 
position vector on the cylinder Cm.n, then ro — rrn,n + u with u = |u| < R. If r is a position vector off the scatterers, 

ore 

00 

Ho{k\r-r,\)^ Y, e'''^'"^--'>^Mku)H,{k\v~TraA) 

where 9 and 9m,n a,re the angles between the x-axis and the vectors u and r — r^^^^ respectively. On the other hand, 
the field inside the cylinder Cm,n is given by 

oo 

^c.(ro) = e™'=^ J2 a,,ne''''Mk'u) (B3) 
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with k' = rick for the index of refraction ric of the cyhnders and the coefficients a^.„ given by 



^U.n rt 



-±-j\_-"E.iU., + r... 



Here e^ = ei cos0 + 63 sin0 is the usual radial vector of polar coordinates. In particular, at/_„ is at most of the order 
of (fci?)~l''l in the limit of thin cylinders. 
It follows that 

/ EUro)G{r\ro)dro = 2TTHe'"'''- V a,,„e*"^'"'"F,(fc|r - r„,,„|) / Mk'u)Mku)udu 

As the integral in the above series is of the order of (fc_R)^l'^l+^, it is then justified to approximate the series by its 
0*'* summand in the limit kR ^ 1, so that 

/ E^{T,)G{v\v,)dv, « i7T^R\'"'''^ao,,Ho{k\r - r™,„|) (B4) 



The value of ao,n may then be recovered by setting u = in Eq.(JB3|. It is -Ecj(ro,n)- This establishes Eq.(JB2| 
everywhere off the scatterers. In particular, the fields are determined by the knowledge of their values £'tj(ro.±i) on 
the cylinders Co^±i alone. 

To determine the values _E(^(ro,±i), let n be either 1 or —1 and, ro.„ be substituted for r in Eq. (|Bl|) . The latter 
equation becomes. 






X] / £'^(ro)i:/o(fc|ro^„-ro|)dro 



(B5) 



where the integral over Co,n has been isolated due to the singularity of its integrand at Tq = ro_„. To approximate 
this particular integral, Eq. (jB3p is used to obtain in the leading order of kR; 

/(■R / 2z / ( kR\ 1 

£;„(ro)fl'o(A:|ro,„ - ro|)dro = 27r£;^(ro,„) / J„{k'u)Ho{ku)udu « 7ri?^£'„(ro^„) f IH f 7 + In 



2 y 2 

where 7 is the Euler constant. All the other integrals in Eq. (jB5|) obey the estimate (jB4p . By taking n successively 
equal to 1 then to —1, the following system is obtained: 



<i>oE^{aei + hes) + ^+E^{-he3) ^ e 

27r5o(fc) 

<^-E^{aei + hes) + <^oEuji-he3) 



i{akx-\-hkz) 



(B6) 



-ihkz 



27rS„{k) 
The functions $01*^^ and <l>~ are 

4.*(<i, ft, k,k,) = ^Yl e'""--f/„(ft|(m T o)e, + teal) 

rrt 

The variants of these functions in Eas. ([TU)) are obtained through the formulas, 

1 °° °° i(a;(fe,+27rm) + |z|fc,,,„) 

J2 e™'=-iIo(fc|r-mei|)= ^ , r^O (B7a) 

m — — 00 m— — 00 ' 
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Relation (jBTa^ can be proved by substituting the plane wave representation of the Hankel function 

into the left side of (jB7a[) and carrying out the integration with respect to K^ followed by an application of the Poisson 
summation formula to evaluate the integral with respect to K^. Relation (jB7bl) is the obtained from (IB7ap in the 
limit r — >■ 0, where for the term m = in the left side of Eq. (jB7ap . the asymptotic expansion of the Hankel function 
for a small argument has to be used. 

When the determinant of the system (|B6p is nonzero, Eq.([2]) has a unique solution. Otherwise, the homogeneous 
Lippmann-Schwinger equation admits nonzero solutions: the bound states. These states can then be arbitrarily 
superposed to obtain the general solution. The various expressions for these bound states in Sections IIIII and IIV Al 
are obtained by applying formula (jB7al) to Eq. (jB2p in the absence of the incident wave, i.e., by omitting the term e'*'''. 



Appendix C: Bound states in the continuums I and II (Complements) 

This section of the Appendix gives some of the technical details omitted in sections [ill Bl and IIV Al First, Eq. ([M)) 
is stablished for the sequence {cm]m=i defined in Ea. (l22]) . In the case under consideration, the sequence {<7z,m} is in 
the order, 

qz-i < qz.i < qz,-2 < qz,2 < •■• 

with equalities occurring when k^ = 0. Since the function f : t ^ (0,oo) i-^- t^^e^^'** is strictly decreasing and 
Cm = fiqz.-m) — f{qz,m), it foUows that Cm > with equality holding only if kx = 0. 

Now suppose that kx ^ ^ and, hence, Cm > 0, m = 1, 2, 3.... To complete the proof of Ea. (l24)) . it suffices to show 
that, 

< e-^'^'' and hm ^^^^ = e-^"'' (CI) 



Cm m-s-oo c„ 



To establish the first of conditions (ICll) . the ratio of Cm+i to Cm is rewritten as. 



Cm+1 a-m-1 ^ a_„_i 6 '^'-■'" 

= 1 a^ : "m = (C2) 

Cm O^—m -L ^_ qz.m 

Next, the following chain of conclusions holds: 

f 1 _ am + l 

\ qz,m + qz,-m-l ^ qz.-m + qz.m+1 Oi-rnOim+1 - _, Q-m-i ^^ , Cm+1 . «-m-l , 2h('o, _™ -o, _„_i 1 

y ' => > 1 =^ — j —^ — < 1 => < < e ^^ • ^-' ' 

I qz,mqz, — m — l ^ qz, — rnqz,m+l Oi—m—lCtm -L ~ ^ ™ C-m OL—m, 



The first of conditions (JC1|) then follows as qz.-m — qz.-m-i < — 27r. The limit in (|C1|) follows from (jC2p and the 
limits, 

hm = e hm = e " 7^ 1 

m— >oo Oi—m, m— >oo Oi—rn 

Second, the formula ([25]) is proved. This is done by a repetitive application of Abel's partial summation formula. 
Let Un be defined for each n = 0, 1, 2, ... by. 



m— 1 



y^ Cm^n sin(27ram) 
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The objective is to show that another expression of — uq is ([23| . By Abel's partial summation formula, 

E, sin(7raTO) sin(7ra(?7i + 1)) 
sm 7ra 

m—l ^ ^ 

oo _. oo 

= cot(7ra) ^ {cm+n - Cm+n+i) sm^{TTam) + X X! i'^m+n - Cm+n+i) sin(27ram) 

m—l m—l 

oo j^ ^ 

= C0t(7ra) ^ (Cm+„ - Cm+n+l) sm^{TTam) + -u„ - -Un+l 
m—l 

Thus, 



w 



Wo 



= {-if+'uN+i + 2cot(™) ^ c„ + 2^(-l)"c„,+„ + {-lf+'c„,+N+i sin^nam), ^N = 1,2,3, 



m—l 



By using the first of conditions (|C1I) it is straightforward that uat+i -^ a.s N —i' oo, and Eg. ([25]) follows. 

Third, the functions k i— > ^„(fc, kx,a) defined in Eg . ([25)) are proved to be monotonically decreasing, i.e., dk^n < 0. 
This derivative reads 



-7^ = ~^^^f TT - E ^ 1 - (-1)" cos(27ram)e-"-«-'=^-^ 1 + n^ 




Now, if i > and n is a positive integer; then e~* (l + i (l + (nTr)"^*^)) < 1. Setting t — mrqz^mk^^ shows that all 
summands are positive and hence 5fe\l/„ < 0, Vti = 1, 2, 3, .... The functions k i— ^ ^n(fc, k^) and k i— >■ ^oo(^, fca:) defined 
in Egs. (l44|) and (l49l) are shown to be decreasing in a similar fashion. 

Lastly, the function k^ i— ^ (/Sqo defined in Eg. (j52p is shown to be monotonically decreasing on (0, tt). To establish this 
fact, note that since 9fc^oo < 0, the implicit function theorem implies that the function k^ i— > kao{kx) is continuously 
differentiable and k^{kx) — —dk^'^ oo{dk'^ oo)~^ ■ Now, 

d'^oo ■^ / 2nm + k^ 27rm + 27r — fc^; \ 

^^ ~ ^1 \^{2T:m + k,f^W ' ^i2TTm + 2TT-k,y-ky ^ 

Hence k^{kx) > 0. By logarithmic differentiation, it follows that, 

V'ooik.) _ fcooC (kl - (2^ - k,)^) - fc, {kl - (2^ - fc,)2) - (27r - fc,) (feg, - fcg) 
Voo(fc.) (fci - A:2) (fc^ - (2^ - A:,)2) 

since < k^ < 2Tr — k^ < koo{kx)- 

Appendix D: Approximations 

In this section of the Appendix, it is outlined how the approximations in Eqs. lfTBl) . ([30]), (|38|) and p7p can be 
obtained. The computations of the wavenumbers in the first three of the latter equations being similar, so only 
Egs. ([5n|) is stablished, and the discussion is finished by proving Eq.(j37]). Suppose first that only one diffraction 
channel is open and (a, k^) is in the set L of Eg. ((25|) . The objective is to approximate the value /c„ G {k^, 2tt — kx) 
such that '^n{kn,kx,a) = for the function ^„ defined in Ea. (P5| . This is done by identifying the leading terms in 
^„ near the wavenumber A;„. 

Since the dielectric cylinders forming the double array are assumed to be thin in comparison to the wavelength 
,i.e., kR ^ 1, the quantity (27r^(,(fc))~^ in the expression of ^n is large. Consequently, the wavenumber fc„ such that 
\l/„(fc„, kx, a) = must be close to the diffraction threshold 2n — kx so that the term q~\i in \l/„ is large enough to 
compensate for the magnitude of {2TrSo{k))^^ . Thus a first approximation for the wavenumber fc„ may be found by 
solving the equation, 

1 1- (-l)"cos(27ra)e-"'^'5'-'-i'=-"' _ 



2nSo{k) 92_i 



25 

which is obtained by keeping only the leading terms in the expression of ^„ near fc„. When (—1)" cos(27ra) — 1, then 
the equation becomes, 

"-=0 (Dl) 



27r(5o(fc) qz-i 

In particular if n is not large, this equation has no roots since as qz.-i becomes smaller, then the second summand 
gets closer to mrk~^ and hence is much smaller than the first summand. This was to be expected since in the case 
(— 1)" cos(27ra) = 1, it was already established that bound states exist only for sufficiently large n. Also, the initial 
integer n at which the wavenumber fc„ exists for (— 1)" cos(27ra) = 1 grows as kR — >■ 0. This makes it impossible 
to provide a good approximation for the exponential term in Eq. ljDip that would allow a perturbative solution of 
the said equation. This complexity disappears when (— 1)" cos(27ra) ^ 1. In this case and for n not too large, the 
approximation g-"-'^<i^,-i>^z r:; 1 is valid and Eq. ljDll) becomes, 

1 1- (-l)"cos(27ra) _ 



2TTSo{k) q^^^i 

Thus, 

fc„ « V(2^ - k.,)^ - 47r2(l - (-!)« cos(27ra))252(2^ _ fc^) (D2) 

The first of Eqs.(j3ni) follows by keeping the first two terms in a series expansion of the right-hand side of Ea. (ID2p in 
powers of (5o(27r — k^)- The distance /i„ in Eqs. (pO)) can then be derived as indicated by system (|27)) . 

Eqs. ([T6|) and (|38|) are obtained through similar treatments of the functions defined in Eqs.(fT5| and (|44|) respectively. 
In particular, for the bound states below the continuum, it appears possible to give only the approximate value of 
the wavenumber fc"*" while the wavenumber k~ eludes the perturbation method due to a complicated dependence of 
its existence condition on the size of the cylinders. Similarly, in the case of two open channels, it turns out to be only 
possible to solve for the wavenumbers k2n+i whose existence is not subject to changes in cylinder sizes. Expressions 
analogous to Eq. (jD2p for the wavenumbers fc+ and k2n+i are. 



y/kl - 167r2,52(fc,) (D3a) 



k2n+i ~ V(27r + A:,)2- 16^2^2(2^ + ^^) (Ogb) 

To establish Eg. ([57)1 . we recall that the point fc2"+i'' is solution to the equation. 



^2n+l ^a; ^ 



V kL+i - (2^ - fcx)^ 2n + 1 
Hence, 

I 



47rfc^(l - 2r^) + 47r2 = (1 - r^)Ak, 



2n+ 1 
where kj,^^^ = {2tt + k^f - Afc for Afc given by Ea. (lD3bt . Thus, 

Ank^il - 2r^) + in^ = n^u{l - r^){27: + k^f , u = R'^{ec - 1)' 
One can then look for a series solution k^ — a^ + aiu + a2U^ + ... This leads to Eq.(|37|). 
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